knitr::opts_chunk$set(echo = TRUE)

cmdstanr::set_cmdstan_path(path = "C:/Users/kueng/.cmdstan/cmdstan-2.35.0")

library(tidyverse)
library(R.utils)
library(wbCorr)
library(readxl)
library(kableExtra)
library(brms)
library(bayesplot)
library(see)
library(beepr)
library(DHARMa)
library(digest)



source(file.path('Functions', 'ReportModels.R'))
source(file.path('Functions', 'PrettyTables.R'))
source(file.path('Functions', 'ReportMeasures.R'))
source(file.path('Functions', 'PrepareData.R'))

report_function_hash <- digest::digest(summarize_brms)
system("shutdown /a")
## [1] 1116
# Set options for analysis
use_mi = FALSE
shutdown = FALSE
report_ordinal = FALSE
report_hurdle = TRUE
do_priorsense = FALSE

options(
  dplyr.print_max = 100, 
  brms.backend = 'cmdstan',
  brms.file_refit = ifelse(use_mi, 'never', 'on_change'),
  brms.file_refit = 'on_change',
  #brms.file_refit = 'always',
  error = function() {
    beepr::beep(sound = 5)
    if (shutdown) {
      system("shutdown /s /t 180")
      quit(save = "no", status = 1)
    }
  }
  , es.use_symbols = TRUE
)


####################### Model parameters #######################

iterations = 2000 # 10'000 per chain to achieve 40'000
warmup = 1000

# NO AR!!!
#corstr = 'ar'
#corstr = 'cosy_couple'
#corstr = 'cosy_couple:user'


################################################################

suffix = as.character(iterations)
df <- openxlsx::read.xlsx(file.path('long.xlsx'))
df_original <- df

df_double <- prepare_data(df, recode_pushing = FALSE, use_mi = use_mi)[[1]]

Constructing scales reshaping data (4field) centering data within and between

summary(df_double$pushing)

Min. 1st Qu. Median Mean 3rd Qu. Max. NA’s 0.0000 0.0000 0.0000 0.3468 0.0000 5.0000 2323

Modelling

# For indistinguishable Dyads
model_rows_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'planPlan',
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb'
  )


model_rows_fixed_ordinal <- c(
  model_rows_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rows_fixed[2:length(model_rows_fixed)]
)

model_rows_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rows_random_ordinal <- c(model_rows_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    "Intercept", 
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Actionplan",
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime"
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
rows_to_pack <- list(
  "Within-Person Effects" = c(2,10),
  "Between-Person Effects" = c(11,17),
  "Random Effects" = c(18, 24), 
  "Additional Parameters" = c(25,25)
  )


rows_to_pack_ordinal <- list(
  "Intercepts" = c(1,6),
  "Within-Person Effects" = c(2+5,10+5),
  "Between-Person Effects" = c(11+5,17+5),
  "Random Effects" = c(18+5, 24+5), 
  "Additional Parameters" = c(25+5,25+6)
  )

HURDLE MODELS

# For indistinguishable Dyads
model_rows_fixed_hu <- c(
    'Intercept', 
    'hu_Intercept',
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'planPlan',
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb',
    
    # HURDLE MODEL
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'hu_persuasion_self_cw', 
    'hu_persuasion_partner_cw', 
    'hu_pressure_self_cw', 
    'hu_pressure_partner_cw', 
    'hu_pushing_self_cw', 
    'hu_pushing_partner_cw', 
    'hu_day', 
    'hu_planPlan',
    'hu_weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'hu_persuasion_self_cb',
    'hu_persuasion_partner_cb',
    'hu_pressure_self_cb',
    'hu_pressure_partner_cb',
    'hu_pushing_self_cb',
    'hu_pushing_partner_cb',
    'hu_weartime_self_cb'
  )

model_rows_fixed_hu_ordinal <- c(
  model_rows_fixed_hu[1:2],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rows_fixed_hu[3:length(model_rows_fixed_hu)]
)


model_rows_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(hu_Intercept)',
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # HURDLE
  'sd(hu_persuasion_self_cw)',
  'sd(hu_persuasion_partner_cw)',
  'sd(hu_pressure_self_cw)',
  'sd(hu_pressure_partner_cw)',
  'sd(hu_pushing_self_cw)',
  'sd(hu_pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rows_random_hu_ordinal <- c(model_rows_random_hu,'disc')
# For indistinguishable Dyads
model_rownames_fixed_hu <- c(
    "Intercept", 
    "Hurdle Intercept",
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Actionplan",
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime",
    
    # HURDLE
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Hu Daily persuasion experienced", 
    "Hu Daily persuasion utilized (partner's view)", # OR partner received
    "Hu Daily pressure experienced", 
    "Hu Daily pressure utilized (partner's view)", 
    "Hu Daily pushing experienced", 
    "Hu Daily pushing utilized (partner's view)", 
    "Hu Day", 
    "Hu Actionplan",
    "Hu Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Hu Mean persuasion experienced", 
    "Hu Mean persuasion utilized (partner's view)", 
    "Hu Mean pressure experienced", 
    "Hu Mean pressure utilized (partner's view)", 
    "Hu Mean pushing experienced", 
    "Hu Mean pushing utilized (partner's view)", 
    "Hu Mean weartime"
  )



model_rownames_fixed_hu_ordinal <- c(
  model_rownames_fixed_hu[1:2],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed_hu[3:length(model_rownames_fixed_hu)]
)



model_rownames_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(Hurdle Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  
  # Hurdle
  "sd(Hu Daily persuasion experienced)", 
  "sd(Hu Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Hu Daily pressure experienced)", 
  "sd(Hu Daily pressure utilized (partner's view))", 
  "sd(Hu Daily pushing experienced)", 
  "sd(Hu Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rownames_random_hu_ordinal <- c(model_rownames_random_hu,'disc')
rows_to_pack_hu <- list(
  "Conditional Within-Person Effects" = c(3,11),
  "Conditional Between-Person Effects" = c(12,18),
  
  "Hurdle Within-Person Effects" = c(19,27),
  "Hurdle Between-Person Effects" = c(28,34),
  
  "Random Effects" = c(35, 48), 
  "Additional Parameters" = c(49,49)
  )

rows_to_pack_hu_ordinal <- list(
  "Intercepts" = c(1,7),
  "Conditional Within-Person Effects" = c(3+5,11+5),
  "Conditional Between-Person Effects" = c(12+5,18+5),
  
  "Hurdle Within-Person Effects" = c(19+5,27+5),
  "Hurdle Between-Person Effects" = c(28+5,34+5),
  
  "Random Effects" = c(35+5, 48+5), 
  "Additional Parameters" = c(49+5,49+6)
  )

Self-Reported MVPA

range(df_double$pa_sub, na.rm = T) 
## [1]   0 720
hist(df_double$pa_sub, breaks = 40) 

hist(log(df_double$pa_sub+00000000001), breaks = 40)

WITHOUT PUSHING

Hurdle Lognormal Model

formula <- bf(
  pa_sub ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    
    plan +
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID),
  
  hu = ~ persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    
    plan +
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID)
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA
  , brms::set_prior("normal(0.5, 2.5)", class = "Intercept", dpar = 'hu') # hurdle part
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_Nopushing", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
pa_sub_digest <- digest::digest(pa_sub)
# summarize with rope range for hurdle part
summary_pa_sub_hurdle <- summarize_brms(
  pa_sub, 
  stats_to_report = c('pd','ROPE'),
  rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) 
## Warning in summarize_brms(pa_sub, stats_to_report = c("pd", "ROPE"), rope_range
## = c(-0.18, : Coefficients were exponentiated. Double check if this was
## intended.
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub$data$pa_sub[pa_sub$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub_continuous <- summarize_brms(
  pa_sub, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = rope_range_continuous,
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) 
## Warning in summarize_brms(pa_sub, stats_to_report = c("CI", "SE", "pd", :
## Coefficients were exponentiated. Double check if this was intended.
# Replace only the ROPE and % in Rope columns for rows with 'Hu'
summary_pa_sub <- summary_pa_sub_continuous

columns_to_replace <- c("ROPE", "inside ROPE")

summary_pa_sub[grepl('Hu', rownames(summary_pa_sub)), columns_to_replace] <- 
  summary_pa_sub_hurdle[grepl('Hu', rownames(summary_pa_sub_hurdle)), columns_to_replace]

# Print the updated dataframe
summary_pa_sub %>%
  print_df(rows_to_pack = rows_to_pack_hu)
exp(Est.) SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 39.09*** 2.57 [34.11, 44.31] 1.000 [0.92, 1.08] 0.000 1.000 1519 2186
Hurdle Intercept 0.27*** 0.04 [ 0.20, 0.36] 1.000 [0.84, 1.20] 0.000 1.002 1267 1972
Conditional Within-Person Effects
Daily persuasion experienced 1.02 0.03 [ 0.97, 1.08] 0.828 [0.92, 1.08] 0.980 1.000 2047 2905
Daily persuasion utilized (partner’s view) 1.02 0.02 [ 0.98, 1.07] 0.866 [0.92, 1.08] 0.992 1.000 2764 3010
Daily pressure experienced 0.91* 0.04 [ 0.82, 1.00] 0.978 [0.92, 1.08] 0.329 1.001 4179 2899
Daily pressure utilized (partner’s view) 0.93 0.04 [ 0.85, 1.02] 0.944 [0.92, 1.08] 0.590 1.000 4011 2995
Daily pushing experienced NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Day 0.96 0.06 [ 0.86, 1.09] 0.723 [0.92, 1.08] 0.722 1.001 5576 2644
Actionplan 1.38*** 0.06 [ 1.27, 1.50] 1.000 [0.92, 1.08] 0.000 1.002 6770 3219
Daily weartime NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.10 0.15 [ 0.83, 1.44] 0.757 [0.92, 1.08] 0.335 1.001 924 1425
Mean persuasion utilized (partner’s view) 1.07 0.15 [ 0.82, 1.41] 0.692 [0.92, 1.08] 0.382 1.003 941 1795
Mean pressure experienced 1.15 0.19 [ 0.84, 1.60] 0.801 [0.92, 1.08] 0.268 1.001 1275 2139
Mean pressure utilized (partner’s view) 0.96 0.16 [ 0.69, 1.33] 0.599 [0.92, 1.08] 0.347 1.002 1199 1938
Mean pushing experienced NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean weartime NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 1.55*** 0.11 [ 1.36, 1.79] 1.000 [0.84, 1.20] 0.000 1.001 2972 2563
Hu Daily persuasion utilized (partner’s view) 1.37*** 0.08 [ 1.23, 1.57] 1.000 [0.84, 1.20] 0.007 1.001 4448 2648
Hu Daily pressure experienced 0.99 0.13 [ 0.74, 1.32] 0.545 [0.84, 1.20] 0.796 1.001 4751 2981
Hu Daily pressure utilized (partner’s view) 1.59** 0.28 [ 1.14, 2.40] 0.997 [0.84, 1.20] 0.048 1.001 4149 2710
Hu Daily pushing experienced NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Day 0.85 0.12 [ 0.64, 1.12] 0.868 [0.84, 1.20] 0.538 1.004 7877 3012
Hu Actionplan 9.74*** 0.88 [ 8.17, 11.74] 1.000 [0.84, 1.20] 0.000 1.005 8410 3039
Hu Daily weartime NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 1.53 0.48 [ 0.80, 2.84] 0.902 [0.84, 1.20] 0.194 1.002 1143 1885
Hu Mean persuasion utilized (partner’s view) 1.58 0.50 [ 0.85, 2.94] 0.923 [0.84, 1.20] 0.169 1.003 1108 1848
Hu Mean pressure experienced 0.34** 0.13 [ 0.16, 0.73] 0.997 [0.84, 1.20] 0.010 1.002 1691 2335
Hu Mean pressure utilized (partner’s view) 0.56 0.22 [ 0.26, 1.20] 0.930 [0.84, 1.20] 0.132 1.002 1612 2154
Hu Mean pushing experienced NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.30 0.04 [0.23, 0.40] NA NA NA 1.003 1023 2027
sd(Hurdle Intercept) 0.73 0.10 [0.55, 0.97] NA NA NA 1.001 1448 2056
sd(Daily persuasion experienced) 0.12 0.02 [0.08, 0.17] NA NA NA 1.001 2526 2989
sd(Daily persuasion utilized (partner’s view)) 0.08 0.02 [0.05, 0.13] NA NA NA 1.000 2678 2238
sd(Daily pressure experienced) 0.07 0.06 [0.00, 0.24] NA NA NA 1.005 1479 1735
sd(Daily pressure utilized (partner’s view)) 0.07 0.05 [0.00, 0.19] NA NA NA 1.002 2054 2288
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion experienced) 0.23 0.08 [0.07, 0.41] NA NA NA 1.002 1302 1082
sd(Hu Daily persuasion utilized (partner’s view)) 0.16 0.09 [0.01, 0.33] NA NA NA 1.002 1146 997
sd(Hu Daily pressure experienced) 0.17 0.15 [0.01, 0.66] NA NA NA 1.003 1846 2077
sd(Hu Daily pressure utilized (partner’s view)) 0.22 0.20 [0.01, 0.86] NA NA NA 1.001 1680 2166
sd(Hu Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma 0.68 0.01 [0.66, 0.70] NA NA NA 1.001 7369 2900
# Plot continuous part of model

variable <- c(
  '(Intercept)',
  'b_persuasion_self_cw',
  'b_persuasion_partner_cw',
  'b_pressure_self_cw',
  'b_pressure_partner_cw'
)


plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(
    pa_sub, 
    parameter = variable, 
    range = rope_range_continuous,
    verbose = F,
    ci = 1
  )
) + theme_bw()

# Hurdle part of the model
variable <- c(
  'b_hu_persuasion_self_cw',
  'b_hu_persuasion_partner_cw',
  'b_hu_pressure_self_cw',
  'b_hu_pressure_partner_cw'
)

plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()

# The rope range for the bernoulli part of the model is -0.18, 0.18
plot(
  bayestestR::rope(pa_sub, parameter = variable, range = c(-0.18, 0.18), ci = 1),
  verbose = FALSE
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.77), b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.79), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.78), b_hu_persuasion_partner_cb and
##   b_hu_persuasion_self_cb (r = 0.73), b_hu_pressure_self_cb and
##   b_hu_persuasion_self_cb (r = 0.76), b_hu_pressure_partner_cb and
##   b_hu_persuasion_partner_cb (r = 0.76). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

Device Based MVPA

Lognormal Model

formula <- bf(
  pa_obj ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    
    plan +
    day + weartime_self_cw + weartime_self_cb +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID)
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 50)", class = "Intercept") 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_nopushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
pa_obj_log_digest <- digest::digest(pa_obj_log)
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summarize_brms(
  pa_obj_log, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = rope_range_log,
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Warning in summarize_brms(pa_obj_log, stats_to_report = c("CI", "SE", "pd", :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 111.65*** 5.82 [101.28, 124.78] 1.000 [0.94, 1.07] 0.000 1.007 717 1643
Within-Person Effects
Daily persuasion experienced 1.03 0.02 [ 1.00, 1.06] 0.973 [0.94, 1.07] 0.990 1.000 2672 3035
Daily persuasion utilized (partner’s view) 1.02 0.02 [ 0.99, 1.05] 0.892 [0.94, 1.07] 0.997 1.000 3291 3157
Daily pressure experienced 0.95 0.03 [ 0.89, 1.02] 0.926 [0.94, 1.07] 0.685 1.001 5011 2874
Daily pressure utilized (partner’s view) 0.98 0.03 [ 0.92, 1.05] 0.691 [0.94, 1.07] 0.930 1.001 4611 3115
Daily pushing experienced NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Day 0.97 0.04 [ 0.90, 1.04] 0.777 [0.94, 1.07] 0.837 1.003 7630 2798
Actionplan 1.09*** 0.02 [ 1.05, 1.14] 1.000 [0.94, 1.07] 0.144 1.000 8010 2892
Daily weartime 1.00*** 0.00 [ 1.00, 1.00] 1.000 [0.94, 1.07] 1.000 1.000 4046 2780
Between-Person Effects
Mean persuasion experienced 1.10 0.14 [ 0.86, 1.42] 0.787 [0.94, 1.07] 0.306 1.005 796 1443
Mean persuasion utilized (partner’s view) 1.02 0.13 [ 0.80, 1.32] 0.573 [0.94, 1.07] 0.392 1.006 784 1435
Mean pressure experienced 0.91 0.13 [ 0.69, 1.18] 0.749 [0.94, 1.07] 0.294 1.004 982 1686
Mean pressure utilized (partner’s view) 1.04 0.14 [ 0.79, 1.35] 0.595 [0.94, 1.07] 0.351 1.008 842 1595
Mean pushing experienced NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean weartime 1.00 0.00 [ 1.00, 1.00] 0.944 [0.94, 1.07] 1.000 1.001 4896 3276
Random Effects
sd(Intercept) 0.29 0.04 [0.23, 0.38] NA NA NA 1.001 876 1677
sd(Daily persuasion experienced) 0.05 0.01 [0.02, 0.08] NA NA NA 1.002 2216 2269
sd(Daily persuasion utilized (partner’s view)) 0.05 0.02 [0.02, 0.09] NA NA NA 1.001 1646 1915
sd(Daily pressure experienced) 0.04 0.04 [0.00, 0.15] NA NA NA 1.000 1899 2148
sd(Daily pressure utilized (partner’s view)) 0.04 0.03 [0.00, 0.11] NA NA NA 1.000 2228 2102
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma 0.58 0.01 [0.56, 0.59] NA NA NA 1.002 5669 3240
plot(
  bayestestR::p_direction(pa_obj_log),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_obj_log, range = rope_range_log, ci = 1)
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.86), b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.84), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.78), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.78), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.86). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

# Nothing significant, no plots

Affect

Gaussian

formula <- bf(
  aff ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +

    plan +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  ,brms::set_prior("normal(0, 20)", class = "Intercept", lb=1, ub=6)
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_nopushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
mood_gauss_digest <- digest::digest(mood_gauss)
summarize_brms(
  mood_gauss, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = F) %>%
  print_df(rows_to_pack = rows_to_pack)
Est. SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 3.66*** 0.11 [ 3.45, 3.86] 1.000 [-0.11, 0.11] 0.000 1.002 349 743
Within-Person Effects
Daily persuasion experienced 0.00 0.02 [-0.04, 0.05] 0.565 [-0.11, 0.11] 1.000 1.001 3204 2834
Daily persuasion utilized (partner’s view) 0.03 0.02 [-0.01, 0.08] 0.912 [-0.11, 0.11] 0.999 1.000 2517 2838
Daily pressure experienced -0.03 0.05 [-0.13, 0.07] 0.721 [-0.11, 0.11] 0.944 1.000 3693 3168
Daily pressure utilized (partner’s view) 0.00 0.05 [-0.10, 0.11] 0.535 [-0.11, 0.11] 0.967 1.001 3267 2612
Daily pushing experienced NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Day 0.25*** 0.06 [ 0.14, 0.37] 1.000 [-0.11, 0.11] 0.009 1.000 6114 2852
Actionplan 0.10** 0.03 [ 0.03, 0.17] 0.999 [-0.11, 0.11] 0.617 1.002 6449 3035
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 0.40 0.25 [-0.10, 0.88] 0.946 [-0.11, 0.11] 0.108 1.003 512 987
Mean persuasion utilized (partner’s view) 0.31 0.25 [-0.18, 0.79] 0.899 [-0.11, 0.11] 0.178 1.003 518 889
Mean pressure experienced -0.32 0.27 [-0.85, 0.20] 0.890 [-0.11, 0.11] 0.169 1.004 575 1114
Mean pressure utilized (partner’s view) -0.26 0.27 [-0.78, 0.27] 0.830 [-0.11, 0.11] 0.219 1.003 581 1179
Mean pushing experienced NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.59 0.08 [0.47, 0.78] NA NA NA 1.008 658 1266
sd(Daily persuasion experienced) 0.05 0.03 [0.00, 0.10] NA NA NA 1.003 906 1079
sd(Daily persuasion utilized (partner’s view)) 0.08 0.03 [0.03, 0.13] NA NA NA 1.001 1496 1601
sd(Daily pressure experienced) 0.06 0.06 [0.00, 0.22] NA NA NA 1.000 1716 1674
sd(Daily pressure utilized (partner’s view)) 0.07 0.06 [0.00, 0.24] NA NA NA 1.002 1396 1589
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma 0.96 0.01 [0.94, 0.98] NA NA NA 1.000 6259 3060
plot(
  bayestestR::p_direction(mood_gauss),
  priors = TRUE
)  + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(mood_gauss, ci = 1)
) + theme_bw()
## Possible multicollinearity between b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.9), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.84), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.85), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.89), b_pressure_partner_cb and
##   b_pressure_self_cb (r = 0.8). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

Reactance

Ordinal

df_double$reactance_ordinal <- factor(df_double$reactance,
                                      levels = 0:5, 
                                      ordered = TRUE)

formula <- bf(
  reactance_ordinal ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +

    plan +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = cumulative() # HURDLE_CUMULATIVE
#  )


#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance_ordinal <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::cumulative(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777
  , file = file.path("models_cache_brms", paste0("reactance_ordinal_nopushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
reactance_ordinal_digest <- digest::digest(reactance_ordinal)
summarize_brms(
  reactance_ordinal, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed_ordinal,
  model_rows_random = model_rows_random_ordinal,
  model_rownames_fixed = model_rownames_fixed_ordinal,
  model_rownames_random = model_rownames_random_ordinal,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_ordinal)
OR SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercepts
Intercept NA NA NA NA NA NA NA NA NA
Intercept[1] 3.58*** 1.04 [ 2.03, 6.38] 1.000 [0.84, 1.20] 0.000 1.000 3499 2918
Intercept[2] 7.49*** 2.25 [ 4.22, 13.84] 1.000 [0.84, 1.20] 0.000 0.999 3448 2962
Intercept[3] 20.32*** 6.43 [ 11.09, 39.16] 1.000 [0.84, 1.20] 0.000 1.000 3542 3144
Intercept[4] 87.96*** 32.35 [ 42.93, 187.99] 1.000 [0.84, 1.20] 0.000 1.001 4026 3073
Intercept[5] 2856.04*** 1828.76 [894.50, 11588.88] 1.000 [0.84, 1.20] 0.000 1.000 5154 3183
Within-Person Effects
Daily persuasion experienced 0.85* 0.07 [ 0.72, 0.99] 0.983 [0.84, 1.20] 0.593 1.001 3783 3010
Daily persuasion utilized (partner’s view) 1.01 0.09 [ 0.84, 1.20] 0.548 [0.84, 1.20] 0.949 1.001 3102 2913
Daily pressure experienced 2.01** 0.35 [ 1.33, 2.82] 0.998 [0.84, 1.20] 0.009 1.001 2513 2764
Daily pressure utilized (partner’s view) 1.15 0.23 [ 0.73, 1.84] 0.761 [0.84, 1.20] 0.516 1.001 2456 2104
Daily pushing experienced NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Day 1.64 0.56 [ 0.84, 3.18] 0.929 [0.84, 1.20] 0.154 1.000 5945 3187
Actionplan 1.03 0.21 [ 0.70, 1.50] 0.551 [0.84, 1.20] 0.618 1.000 4672 2843
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.08 0.51 [ 0.42, 2.80] 0.572 [0.84, 1.20] 0.302 1.004 1403 2385
Mean persuasion utilized (partner’s view) 0.78 0.41 [ 0.30, 2.20] 0.677 [0.84, 1.20] 0.238 1.005 1606 2410
Mean pressure experienced 3.32* 1.70 [ 1.19, 9.62] 0.985 [0.84, 1.20] 0.019 1.004 1692 2235
Mean pressure utilized (partner’s view) 1.23 0.68 [ 0.40, 3.42] 0.647 [0.84, 1.20] 0.236 1.004 1685 2580
Mean pushing experienced NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.78 0.19 [0.47, 1.22] NA NA NA 1.002 1599 2400
sd(Daily persuasion experienced) 0.17 0.13 [0.01, 0.43] NA NA NA 1.006 638 1049
sd(Daily persuasion utilized (partner’s view)) 0.20 0.12 [0.02, 0.48] NA NA NA 1.004 1126 1063
sd(Daily pressure experienced) 0.47 0.22 [0.07, 1.01] NA NA NA 1.003 1001 847
sd(Daily pressure utilized (partner’s view)) 0.29 0.29 [0.01, 1.41] NA NA NA 1.001 1145 2311
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA
disc 1.00 0.00 [1.00, 1.00] NA NA NA NA NA NA
plot(
  bayestestR::p_direction(reactance_ordinal),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()

plot(
  bayestestR::rope(reactance_ordinal, range = c(-0.18, 0.18), ci = 1)
) + theme_bw()
## Possible multicollinearity between b_Intercept[4] and b_Intercept[2] (r
##   = 0.8), b_Intercept[4] and b_Intercept[3] (r = 0.86), b_pressure_self_cb
##   and b_persuasion_self_cb (r = 0.8), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.84). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

conditional_spaghetti(
  reactance_ordinal, 
  effects = c('persuasion_self_cw', 'pressure_self_cw')
  , group_var = 'coupleID'
  #, n_groups = 15
  , plot_full_range = T
)

\(persuasion_self_cw <img src="03_SensitivityPushingSeparate_files/figure-html/report_reactance_ordinal-3.png" width="2400" />\)pressure_self_cw

Binary

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}



df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}


formula <- bf(
  is_reactance ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +

    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +

    plan +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw | coupleID)
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 10)", class = "Intercept", lb=0, ub=5) 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_nopushing_", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
is_reactance_digest <- digest::digest(is_reactance)
summarize_brms(
  is_reactance, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
OR SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 0.31*** 0.10 [0.16, 0.57] 1.000 [0.83, 1.20] 0.000 1.001 4202 3377
Within-Person Effects
Daily persuasion experienced 0.85 0.08 [0.71, 1.02] 0.957 [0.83, 1.20] 0.605 1.000 3643 3198
Daily persuasion utilized (partner’s view) 1.10 0.14 [0.86, 1.46] 0.793 [0.83, 1.20] 0.719 1.000 3319 3084
Daily pressure experienced 2.09* 0.53 [1.25, 4.16] 0.994 [0.83, 1.20] 0.020 1.001 2178 1389
Daily pressure utilized (partner’s view) 1.33 0.50 [0.63, 4.00] 0.786 [0.83, 1.20] 0.287 1.003 1911 1442
Daily pushing experienced NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Day 1.84 0.68 [0.90, 3.72] 0.957 [0.83, 1.20] 0.107 1.001 6332 2915
Actionplan 1.07 0.25 [0.67, 1.74] 0.622 [0.83, 1.20] 0.542 0.999 5977 3011
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.91 1.12 [0.63, 6.06] 0.871 [0.83, 1.20] 0.142 1.002 1488 1842
Mean persuasion utilized (partner’s view) 1.07 0.66 [0.33, 3.87] 0.539 [0.83, 1.20] 0.226 1.001 1732 2699
Mean pressure experienced 17.22*** 15.48 [3.14, 113.96] 1.000 [0.83, 1.20] 0.001 1.001 2115 2134
Mean pressure utilized (partner’s view) 1.15 1.09 [0.17, 7.83] 0.556 [0.83, 1.20] 0.146 1.002 2328 3018
Mean pushing experienced NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 1.21 0.24 [0.81, 1.76] NA NA NA 1.004 1669 2622
sd(Daily persuasion experienced) 0.19 0.13 [0.01, 0.50] NA NA NA 1.008 908 1526
sd(Daily persuasion utilized (partner’s view)) 0.43 0.17 [0.12, 0.84] NA NA NA 1.002 1349 988
sd(Daily pressure experienced) 0.77 0.53 [0.05, 2.08] NA NA NA 1.004 775 1327
sd(Daily pressure utilized (partner’s view)) 0.74 0.66 [0.03, 2.59] NA NA NA 1.006 1088 1559
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA
plot(
  bayestestR::p_direction(is_reactance),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(is_reactance, ci = 1)
) + theme_bw()

Report All Models

process_model_component <- function(obj) {
  # Convert to string, modify, and evaluate
  name <- deparse(substitute(obj))
  if (report_hurdle) name <- paste0(name, '_hu')
  if (report_ordinal) name <- paste0(name, '_ordinal')
  return(get(name, envir = parent.frame()))
}



all_models <- report_side_by_side(
  pa_sub,
  pa_obj_log,
  mood_gauss,
  reactance_ordinal,
  is_reactance,
  
  stats_to_report = c('CI', 'pd'),
  
  model_rows_random = process_model_component(model_rows_random),
  model_rows_fixed = process_model_component(model_rows_fixed),
  model_rownames_random = process_model_component(model_rownames_random),
  model_rownames_fixed = process_model_component(model_rownames_fixed)
) 

[1] “pa_sub”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “pa_obj_log”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “mood_gauss” [1] “reactance_ordinal” [1] “is_reactance”

# pretty printing

summary_all_models <- all_models %>%
  print_df(rows_to_pack = process_model_component(rows_to_pack)) %>%
  add_header_above(
    c(" ", "Subjective MVPA Hurdle Lognormal" = (length(all_models) / 5),  
      "Device-Based MVPA Log (Gaussian)" = (length(all_models) / 5), 
      "Mood Gaussian" = (length(all_models) / 5),
      "Reactance Ordinal" = (length(all_models) / 5),
      "Reactance Dichotome" = (length(all_models) / 5))
  )


export_xlsx(
  summary_all_models, 
  rows_to_pack = process_model_component(rows_to_pack),
  file.path("Output", "AllModelsNoPushing.xlsx"), 
  merge_option = 'header', 
  simplify_2nd_row = TRUE,
  line_above_rows = c(1,2),
  line_below_rows = c(-1)
)
## 
## Attaching package: 'rvest'
## The following object is masked from 'package:readr':
## 
##     guess_encoding
print(summary_all_models)
Subjective MVPA Hurdle Lognormal
Device-Based MVPA Log (Gaussian)
Mood Gaussian
Reactance Ordinal
Reactance Dichotome
exp(Est.) pa_sub 95% CI pa_sub pd pa_sub exp(Est.) pa_obj_log 95% CI pa_obj_log pd pa_obj_log Est. mood_gauss 95% CI mood_gauss pd mood_gauss OR reactance_ordinal 95% CI reactance_ordinal pd reactance_ordinal OR is_reactance 95% CI is_reactance pd is_reactance
Intercept 39.09*** [34.11, 44.31] 1.000 111.65*** [101.28, 124.78] 1.000 3.66*** [ 3.45, 3.86] 1.000 NA NA NA 0.31*** [0.16, 0.57] 1.000
Hurdle Intercept 0.27*** [ 0.20, 0.36] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced 1.02 [ 0.97, 1.08] 0.828 1.03 [ 1.00, 1.06] 0.973 0.00 [-0.04, 0.05] 0.565 0.85* [ 0.72, 0.99] 0.983 0.85 [0.71, 1.02] 0.957
Daily persuasion utilized (partner’s view) 1.02 [ 0.98, 1.07] 0.866 1.02 [ 0.99, 1.05] 0.892 0.03 [-0.01, 0.08] 0.912 1.01 [ 0.84, 1.20] 0.548 1.10 [0.86, 1.46] 0.793
Daily pressure experienced 0.91* [ 0.82, 1.00] 0.978 0.95 [ 0.89, 1.02] 0.926 -0.03 [-0.13, 0.07] 0.721 2.01** [ 1.33, 2.82] 0.998 2.09* [1.25, 4.16] 0.994
Daily pressure utilized (partner’s view) 0.93 [ 0.85, 1.02] 0.944 0.98 [ 0.92, 1.05] 0.691 0.00 [-0.10, 0.11] 0.535 1.15 [ 0.73, 1.84] 0.761 1.33 [0.63, 4.00] 0.786
Daily pushing experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Day 0.96 [ 0.86, 1.09] 0.723 0.97 [ 0.90, 1.04] 0.777 0.25*** [ 0.14, 0.37] 1.000 1.64 [ 0.84, 3.18] 0.929 1.84 [0.90, 3.72] 0.957
Actionplan 1.38*** [ 1.27, 1.50] 1.000 1.09*** [ 1.05, 1.14] 1.000 0.10** [ 0.03, 0.17] 0.999 1.03 [ 0.70, 1.50] 0.551 1.07 [0.67, 1.74] 0.622
Daily weartime NA NA NA 1.00*** [ 1.00, 1.00] 1.000 NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.10 [ 0.83, 1.44] 0.757 1.10 [ 0.86, 1.42] 0.787 0.40 [-0.10, 0.88] 0.946 1.08 [ 0.42, 2.80] 0.572 1.91 [0.63, 6.06] 0.871
Mean persuasion utilized (partner’s view) 1.07 [ 0.82, 1.41] 0.692 1.02 [ 0.80, 1.32] 0.573 0.31 [-0.18, 0.79] 0.899 0.78 [ 0.30, 2.20] 0.677 1.07 [0.33, 3.87] 0.539
Mean pressure experienced 1.15 [ 0.84, 1.60] 0.801 0.91 [ 0.69, 1.18] 0.749 -0.32 [-0.85, 0.20] 0.890 3.32* [ 1.19, 9.62] 0.985 17.22*** [3.14, 113.96] 1.000
Mean pressure utilized (partner’s view) 0.96 [ 0.69, 1.33] 0.599 1.04 [ 0.79, 1.35] 0.595 -0.26 [-0.78, 0.27] 0.830 1.23 [ 0.40, 3.42] 0.647 1.15 [0.17, 7.83] 0.556
Mean pushing experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean weartime NA NA NA 1.00 [ 1.00, 1.00] 0.944 NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 1.55*** [ 1.36, 1.79] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) 1.37*** [ 1.23, 1.57] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced 0.99 [ 0.74, 1.32] 0.545 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) 1.59** [ 1.14, 2.40] 0.997 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Day 0.85 [ 0.64, 1.12] 0.868 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Actionplan 9.74*** [ 8.17, 11.74] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 1.53 [ 0.80, 2.84] 0.902 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) 1.58 [ 0.85, 2.94] 0.923 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced 0.34** [ 0.16, 0.73] 0.997 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) 0.56 [ 0.26, 1.20] 0.930 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.30 [0.23, 0.40] NA 0.29 [0.23, 0.38] NA 0.59 [0.47, 0.78] NA 0.78 [0.47, 1.22] NA 1.21 [0.81, 1.76] NA
sd(Hurdle Intercept) 0.73 [0.55, 0.97] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily persuasion experienced) 0.12 [0.08, 0.17] NA 0.05 [0.02, 0.08] NA 0.05 [0.00, 0.10] NA 0.17 [0.01, 0.43] NA 0.19 [0.01, 0.50] NA
sd(Daily persuasion utilized (partner’s view)) 0.08 [0.05, 0.13] NA 0.05 [0.02, 0.09] NA 0.08 [0.03, 0.13] NA 0.20 [0.02, 0.48] NA 0.43 [0.12, 0.84] NA
sd(Daily pressure experienced) 0.07 [0.00, 0.24] NA 0.04 [0.00, 0.15] NA 0.06 [0.00, 0.22] NA 0.47 [0.07, 1.01] NA 0.77 [0.05, 2.08] NA
sd(Daily pressure utilized (partner’s view)) 0.07 [0.00, 0.19] NA 0.04 [0.00, 0.11] NA 0.07 [0.00, 0.24] NA 0.29 [0.01, 1.41] NA 0.74 [0.03, 2.59] NA
sd(Daily pushing experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion experienced) 0.23 [0.07, 0.41] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.16 [0.01, 0.33] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) 0.17 [0.01, 0.66] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) 0.22 [0.01, 0.86] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma 0.68 [0.66, 0.70] NA 0.58 [0.56, 0.59] NA 0.96 [0.94, 0.98] NA NA NA NA NA NA NA

ONLY PUSHING

Hurdle Lognormal Model

formula <- bf(
  pa_sub ~ 
    pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day + 
    
    # Random effects
    (1 + pushing_self_cw + pushing_partner_cw | coupleID),
  
  hu = ~ pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day + 
    
    # Random effects
    (1 + pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA
  , brms::set_prior("normal(0.5, 2.5)", class = "Intercept", dpar = 'hu') # hurdle part
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_onlypushing", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
pa_sub_digest <- digest::digest(pa_sub)
# summarize with rope range for hurdle part
summary_pa_sub_hurdle <- summarize_brms(
  pa_sub, 
  stats_to_report = c('pd','ROPE'),
  rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) 
## Warning in summarize_brms(pa_sub, stats_to_report = c("pd", "ROPE"), rope_range
## = c(-0.18, : Coefficients were exponentiated. Double check if this was
## intended.
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub$data$pa_sub[pa_sub$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub_continuous <- summarize_brms(
  pa_sub, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = rope_range_continuous,
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) 
## Warning in summarize_brms(pa_sub, stats_to_report = c("CI", "SE", "pd", :
## Coefficients were exponentiated. Double check if this was intended.
# Replace only the ROPE and % in Rope columns for rows with 'Hu'
summary_pa_sub <- summary_pa_sub_continuous

columns_to_replace <- c("ROPE", "inside ROPE")

summary_pa_sub[grepl('Hu', rownames(summary_pa_sub)), columns_to_replace] <- 
  summary_pa_sub_hurdle[grepl('Hu', rownames(summary_pa_sub_hurdle)), columns_to_replace]

# Print the updated dataframe
summary_pa_sub %>%
  print_df(rows_to_pack = rows_to_pack_hu)
exp(Est.) SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 57.28*** 4.53 [49.10, 66.97] 1.000 [0.93, 1.08] 0.000 1.002 1327 2040
Hurdle Intercept 3.79*** 0.98 [ 2.27, 6.26] 1.000 [0.84, 1.20] 0.000 1.005 1088 2463
Conditional Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pushing experienced 0.97 0.03 [ 0.91, 1.03] 0.856 [0.93, 1.08] 0.931 1.002 4427 2662
Daily pushing utilized (partner’s view) 0.97 0.03 [ 0.90, 1.03] 0.847 [0.93, 1.08] 0.907 1.001 3679 2873
Day 0.91 0.07 [ 0.79, 1.06] 0.894 [0.93, 1.08] 0.398 1.000 8744 2600
Actionplan NA NA NA NA NA NA NA NA NA
Daily weartime NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pushing experienced 1.32* 0.18 [ 1.00, 1.71] 0.977 [0.93, 1.08] 0.060 1.005 1217 2176
Mean pushing utilized (partner’s view) 1.18 0.15 [ 0.90, 1.53] 0.883 [0.93, 1.08] 0.222 1.005 1265 2076
Mean weartime NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced 1.55** 0.23 [ 1.17, 2.24] 0.998 [0.84, 1.20] 0.035 1.002 3127 2531
Hu Daily pushing utilized (partner’s view) 1.84*** 0.30 [ 1.38, 2.68] 1.000 [0.84, 1.20] 0.002 1.001 2766 2989
Hu Day 0.61 0.15 [ 0.37, 1.00] 0.974 [0.84, 1.20] 0.100 1.003 7632 2768
Hu Actionplan NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced 0.46** 0.14 [ 0.26, 0.81] 0.997 [0.84, 1.20] 0.018 1.001 2675 2824
Hu Mean pushing utilized (partner’s view) 0.62 0.18 [ 0.33, 1.08] 0.952 [0.84, 1.20] 0.136 1.001 2568 2699
Hu Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.37 0.06 [0.29, 0.50] NA NA NA 1.002 1283 2370
sd(Hurdle Intercept) 1.25 0.19 [0.93, 1.74] NA NA NA 1.002 1448 2426
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.06 0.04 [0.00, 0.16] NA NA NA 1.001 870 1486
sd(Daily pushing utilized (partner’s view)) 0.09 0.04 [0.01, 0.18] NA NA NA 1.004 977 1142
sd(Hu Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) 0.46 0.24 [0.07, 1.02] NA NA NA 1.000 947 1689
sd(Hu Daily pushing utilized (partner’s view)) 0.51 0.24 [0.10, 1.13] NA NA NA 1.001 1088 1582
Additional Parameters
sigma 0.67 0.02 [0.64, 0.71] NA NA NA 1.002 6374 2773
# Plot continuous part of model

variable <- c(
  '(Intercept)',
  'b_pushing_self_cw',
  'b_pushing_partner_cw'
)


plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()

plot(
  bayestestR::rope(
    pa_sub, 
    parameter = variable, 
    range = rope_range_continuous,
    verbose = F,
    ci = 1
  )
) + theme_bw()

# Hurdle part of the model
variable <- c(
  'b_hu_pushing_self_cw',
  'b_hu_pushing_partner_cw'
)

plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()

# The rope range for the bernoulli part of the model is -0.18, 0.18
plot(
  bayestestR::rope(pa_sub, parameter = variable, range = c(-0.18, 0.18), ci = 1),
  verbose = FALSE
) + theme_bw()
## Possible multicollinearity between b_pushing_partner_cb and
##   b_pushing_self_cb (r = 0.77). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

Device Based MVPA

Lognormal Model

formula <- bf(
  pa_obj ~ 
   
    pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day + weartime_self_cw + weartime_self_cb +
    
    # Random effects
    (pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 50)", class = "Intercept") 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_onlypushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
pa_obj_log_digest <- digest::digest(pa_obj_log)
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summarize_brms(
  pa_obj_log, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = rope_range_log,
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Warning in summarize_brms(pa_obj_log, stats_to_report = c("CI", "SE", "pd", :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 124.24*** 7.52 [109.69, 139.58] 1.000 [0.94, 1.06] 0.000 1.002 967 1634
Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pushing experienced 1.03 0.03 [ 0.97, 1.08] 0.807 [0.94, 1.06] 0.914 1.001 3221 2733
Daily pushing utilized (partner’s view) 1.02 0.02 [ 0.97, 1.06] 0.797 [0.94, 1.06] 0.979 1.000 5065 2748
Day 0.97 0.05 [ 0.87, 1.08] 0.718 [0.94, 1.06] 0.654 1.000 7350 3336
Actionplan NA NA NA NA NA NA NA NA NA
Daily weartime 1.00** 0.00 [ 1.00, 1.00] 0.998 [0.94, 1.06] 1.000 1.001 4680 2898
Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pushing experienced 0.97 0.06 [ 0.86, 1.10] 0.680 [0.94, 1.06] 0.635 1.001 2223 2645
Mean pushing utilized (partner’s view) 1.04 0.07 [ 0.92, 1.18] 0.739 [0.94, 1.06] 0.578 1.002 2058 2692
Mean weartime 1.00 0.00 [ 1.00, 1.00] 0.812 [0.94, 1.06] 1.000 1.000 3879 3246
Random Effects
sd(Intercept) 0.31 0.04 [0.24, 0.42] NA NA NA 1.005 890 1790
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.09 0.04 [0.02, 0.17] NA NA NA 1.004 887 821
sd(Daily pushing utilized (partner’s view)) 0.03 0.03 [0.00, 0.10] NA NA NA 1.001 1398 1731
Additional Parameters
sigma 0.54 0.01 [0.52, 0.57] NA NA NA 1.001 6040 3037
plot(
  bayestestR::p_direction(pa_obj_log),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_obj_log, range = rope_range_log, ci = 1)
) + theme_bw()

# Nothing significant, no plots

Affect

Gaussian

formula <- bf(
  aff ~ 
    pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day +
    
    # Random effects
    (pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  ,brms::set_prior("normal(0, 20)", class = "Intercept", lb=1, ub=6)
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_onlypushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
mood_gauss_digest <- digest::digest(mood_gauss)
summarize_brms(
  mood_gauss, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = F) %>%
  print_df(rows_to_pack = rows_to_pack)
Est. SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 3.78*** 0.12 [ 3.56, 4.00] 1.000 [-0.11, 0.11] 0.000 1.013 560 1111
Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pushing experienced 0.02 0.04 [-0.05, 0.08] 0.670 [-0.11, 0.11] 0.995 1.002 2823 2591
Daily pushing utilized (partner’s view) 0.09* 0.04 [ 0.01, 0.17] 0.984 [-0.11, 0.11] 0.715 1.002 1597 2207
Day 0.28*** 0.09 [ 0.11, 0.46] 1.000 [-0.11, 0.11] 0.023 1.001 3871 3088
Actionplan NA NA NA NA NA NA NA NA NA
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pushing experienced 0.12 0.11 [-0.08, 0.32] 0.880 [-0.11, 0.11] 0.433 1.003 1269 1954
Mean pushing utilized (partner’s view) 0.08 0.11 [-0.12, 0.30] 0.786 [-0.11, 0.11] 0.554 1.002 1216 1958
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.61 0.08 [0.48, 0.79] NA NA NA 1.011 648 1448
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.07 0.04 [0.00, 0.16] NA NA NA 1.000 1409 1126
sd(Daily pushing utilized (partner’s view)) 0.14 0.04 [0.06, 0.24] NA NA NA 1.002 1329 828
Additional Parameters
sigma 0.91 0.02 [0.87, 0.94] NA NA NA 1.001 3887 2536
plot(
  bayestestR::p_direction(mood_gauss),
  priors = TRUE
)  + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(mood_gauss, ci = 1)
) + theme_bw()

conditional_spaghetti(
  mood_gauss, 
  effects = c('pushing_partner_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

$pushing_partner_cw

Reactance

Ordinal

df_double$reactance_ordinal <- factor(df_double$reactance,
                                      levels = 0:5, 
                                      ordered = TRUE)

formula <- bf(
  reactance_ordinal ~ 
    pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day +
    
    # Random effects
    (pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = cumulative() # HURDLE_CUMULATIVE
#  )


#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance_ordinal <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::cumulative(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777
  , file = file.path("models_cache_brms", paste0("reactance_ordinal_onlypushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
reactance_ordinal_digest <- digest::digest(reactance_ordinal)
summarize_brms(
  reactance_ordinal, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed_ordinal,
  model_rows_random = model_rows_random_ordinal,
  model_rownames_fixed = model_rownames_fixed_ordinal,
  model_rownames_random = model_rownames_random_ordinal,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_ordinal)
OR SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercepts
Intercept NA NA NA NA NA NA NA NA NA
Intercept[1] 5.44*** 1.57 [ 3.06, 10.06] 1.000 [0.84, 1.20] 0.000 1.001 3521 2882
Intercept[2] 9.89*** 3.14 [ 5.45, 18.74] 1.000 [0.84, 1.20] 0.000 1.002 3598 2951
Intercept[3] 23.10*** 7.89 [ 12.38, 46.79] 1.000 [0.84, 1.20] 0.000 1.001 3933 2899
Intercept[4] 81.88*** 35.09 [ 38.53, 199.63] 1.000 [0.84, 1.20] 0.000 1.000 4569 3054
Intercept[5] 771.25*** 625.85 [202.52, 4855.18] 1.000 [0.84, 1.20] 0.000 1.001 5246 2545
Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pushing experienced 1.31* 0.14 [ 1.04, 1.62] 0.987 [0.84, 1.20] 0.214 1.001 3641 2884
Daily pushing utilized (partner’s view) 0.96 0.13 [ 0.75, 1.25] 0.615 [0.84, 1.20] 0.810 1.002 4044 2794
Day 1.60 0.73 [ 0.68, 3.95] 0.854 [0.84, 1.20] 0.182 1.002 4822 3245
Actionplan NA NA NA NA NA NA NA NA NA
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pushing experienced 2.32** 0.65 [ 1.35, 3.99] 1.000 [0.84, 1.20] 0.005 1.001 3394 2987
Mean pushing utilized (partner’s view) 1.25 0.35 [ 0.70, 2.16] 0.783 [0.84, 1.20] 0.365 1.001 3738 2562
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.85 0.21 [0.49, 1.34] NA NA NA 1.001 1569 2204
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.17 0.14 [0.01, 0.50] NA NA NA 1.004 1241 1651
sd(Daily pushing utilized (partner’s view)) 0.16 0.14 [0.01, 0.57] NA NA NA 1.003 1626 1714
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA
disc 1.00 0.00 [1.00, 1.00] NA NA NA NA NA NA
plot(
  bayestestR::p_direction(reactance_ordinal),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()

plot(
  bayestestR::rope(reactance_ordinal, range = c(-0.18, 0.18), ci = 1)
) + theme_bw()
## Possible multicollinearity between b_Intercept[4] and b_Intercept[2] (r
##   = 0.75), b_Intercept[4] and b_Intercept[3] (r = 0.82). This might lead
##   to inappropriate results. See 'Details' in '?rope'.

Binary

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}



df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}


formula <- bf(
  is_reactance ~ 
    pushing_self_cw + pushing_partner_cw +
    pushing_self_cb + pushing_partner_cb +
    
    #plan +
    day +
    
    # Random effects
    (pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 10)", class = "Intercept", lb=0, ub=5) 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_onlypushing", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
is_reactance_digest <- digest::digest(is_reactance)
summarize_brms(
  is_reactance, 
  stats_to_report = c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS'),
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
OR SE 95% CI pd ROPE inside ROPE Rhat Bulk_ESS Tail_ESS
Intercept 0.40** 0.11 [0.23, 0.68] 1.000 [0.83, 1.20] 0.002 1.000 3095 2830
Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Daily pushing experienced 1.37* 0.18 [1.06, 1.87] 0.993 [0.83, 1.20] 0.154 1.000 3238 2410
Daily pushing utilized (partner’s view) 1.16 0.26 [0.80, 2.12] 0.769 [0.83, 1.20] 0.509 1.000 2310 2170
Day 1.78 0.89 [0.69, 4.81] 0.874 [0.83, 1.20] 0.160 1.000 6984 2927
Actionplan NA NA NA NA NA NA NA NA NA
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Mean pushing experienced 5.39*** 2.58 [2.38, 14.94] 1.000 [0.83, 1.20] 0.000 1.001 2293 2628
Mean pushing utilized (partner’s view) 2.37 1.15 [0.97, 6.82] 0.972 [0.83, 1.20] 0.057 1.001 2556 2640
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 1.52 0.31 [1.03, 2.24] NA NA NA 1.001 1431 2494
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.18 0.15 [0.01, 0.56] NA NA NA 1.003 1256 1876
sd(Daily pushing utilized (partner’s view)) 0.42 0.27 [0.03, 1.15] NA NA NA 1.002 1292 1879
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA
plot(
  bayestestR::p_direction(is_reactance),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(is_reactance, ci = 1)
) + theme_bw()

conditional_spaghetti(
  is_reactance, 
  effects = c('pushing_self_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

$pushing_self_cw

Report All Models

process_model_component <- function(obj) {
  # Convert to string, modify, and evaluate
  name <- deparse(substitute(obj))
  if (report_hurdle) name <- paste0(name, '_hu')
  if (report_ordinal) name <- paste0(name, '_ordinal')
  return(get(name, envir = parent.frame()))
}



all_models <- report_side_by_side(
  pa_sub,
  pa_obj_log,
  mood_gauss,
  reactance_ordinal,
  is_reactance,
  
  stats_to_report = c('CI', 'pd'),
  
  model_rows_random = process_model_component(model_rows_random),
  model_rows_fixed = process_model_component(model_rows_fixed),
  model_rownames_random = process_model_component(model_rownames_random),
  model_rownames_fixed = process_model_component(model_rownames_fixed)
) 

[1] “pa_sub”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “pa_obj_log”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “mood_gauss” [1] “reactance_ordinal” [1] “is_reactance”

# pretty printing

summary_all_models <- all_models %>%
  print_df(rows_to_pack = process_model_component(rows_to_pack)) %>%
  add_header_above(
    c(" ", "Subjective MVPA Hurdle Lognormal" = (length(all_models) / 5),  
      "Device-Based MVPA Log (Gaussian)" = (length(all_models) / 5), 
      "Mood Gaussian" = (length(all_models) / 5),
      "Reactance Ordinal" = (length(all_models) / 5),
      "Reactance Dichotome" = (length(all_models) / 5))
  )


export_xlsx(
  summary_all_models, 
  rows_to_pack = process_model_component(rows_to_pack),
  file.path("Output", "AllModelsOnlyPushing.xlsx"), 
  merge_option = 'header', 
  simplify_2nd_row = TRUE,
  line_above_rows = c(1,2),
  line_below_rows = c(-1)
)

  
print(summary_all_models)
Subjective MVPA Hurdle Lognormal
Device-Based MVPA Log (Gaussian)
Mood Gaussian
Reactance Ordinal
Reactance Dichotome
exp(Est.) pa_sub 95% CI pa_sub pd pa_sub exp(Est.) pa_obj_log 95% CI pa_obj_log pd pa_obj_log Est. mood_gauss 95% CI mood_gauss pd mood_gauss OR reactance_ordinal 95% CI reactance_ordinal pd reactance_ordinal OR is_reactance 95% CI is_reactance pd is_reactance
Intercept 57.28*** [49.10, 66.97] 1.000 124.24*** [109.69, 139.58] 1.000 3.78*** [ 3.56, 4.00] 1.000 NA NA NA 0.40** [0.23, 0.68] 1.000
Hurdle Intercept 3.79*** [ 2.27, 6.26] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily pressure experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily pushing experienced 0.97 [ 0.91, 1.03] 0.856 1.03 [ 0.97, 1.08] 0.807 0.02 [-0.05, 0.08] 0.670 1.31* [ 1.04, 1.62] 0.987 1.37* [1.06, 1.87] 0.993
Daily pushing utilized (partner’s view) 0.97 [ 0.90, 1.03] 0.847 1.02 [ 0.97, 1.06] 0.797 0.09* [ 0.01, 0.17] 0.984 0.96 [ 0.75, 1.25] 0.615 1.16 [0.80, 2.12] 0.769
Day 0.91 [ 0.79, 1.06] 0.894 0.97 [ 0.87, 1.08] 0.718 0.28*** [ 0.11, 0.46] 1.000 1.60 [ 0.68, 3.95] 0.854 1.78 [0.69, 4.81] 0.874
Actionplan NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Daily weartime NA NA NA 1.00** [ 1.00, 1.00] 0.998 NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean pressure experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Mean pushing experienced 1.32* [ 1.00, 1.71] 0.977 0.97 [ 0.86, 1.10] 0.680 0.12 [-0.08, 0.32] 0.880 2.32** [ 1.35, 3.99] 1.000 5.39*** [2.38, 14.94] 1.000
Mean pushing utilized (partner’s view) 1.18 [ 0.90, 1.53] 0.883 1.04 [ 0.92, 1.18] 0.739 0.08 [-0.12, 0.30] 0.786 1.25 [ 0.70, 2.16] 0.783 2.37 [0.97, 6.82] 0.972
Mean weartime NA NA NA 1.00 [ 1.00, 1.00] 0.812 NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced 1.55** [ 1.17, 2.24] 0.998 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) 1.84*** [ 1.38, 2.68] 1.000 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Day 0.61 [ 0.37, 1.00] 0.974 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Actionplan NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced 0.46** [ 0.26, 0.81] 0.997 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) 0.62 [ 0.33, 1.08] 0.952 NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.37 [0.29, 0.50] NA 0.31 [0.24, 0.42] NA 0.61 [0.48, 0.79] NA 0.85 [0.49, 1.34] NA 1.52 [1.03, 2.24] NA
sd(Hurdle Intercept) 1.25 [0.93, 1.74] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily persuasion experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily pressure experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily pushing experienced) 0.06 [0.00, 0.16] NA 0.09 [0.02, 0.17] NA 0.07 [0.00, 0.16] NA 0.17 [0.01, 0.50] NA 0.18 [0.01, 0.56] NA
sd(Daily pushing utilized (partner’s view)) 0.09 [0.01, 0.18] NA 0.03 [0.00, 0.10] NA 0.14 [0.06, 0.24] NA 0.16 [0.01, 0.57] NA 0.42 [0.03, 1.15] NA
sd(Hu Daily persuasion experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) 0.46 [0.07, 1.02] NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) 0.51 [0.10, 1.13] NA NA NA NA NA NA NA NA NA NA NA NA NA
Additional Parameters
sigma 0.67 [0.64, 0.71] NA 0.54 [0.52, 0.57] NA 0.91 [0.87, 0.94] NA NA NA NA NA NA NA
report::report_system()

Analyses were conducted using the R Statistical language (version 4.4.1; R Core Team, 2024) on Windows 11 x64 (build 22635)

report::report_packages()
  • beepr (version 2.0; Bååth R, 2024)
  • R.methodsS3 (version 1.8.2; Bengtsson H, 2003)
  • R.oo (version 1.27.0; Bengtsson H, 2003)
  • R.utils (version 2.12.3; Bengtsson H, 2023)
  • brms (version 2.22.0; Bürkner P, 2017)
  • digest (version 0.6.37; Eddelbuettel D, 2024)
  • Rcpp (version 1.0.13.1; Eddelbuettel D et al., 2024)
  • bayesplot (version 1.11.1; Gabry J, Mahr T, 2024)
  • lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  • DHARMa (version 0.4.7; Hartig F, 2024)
  • wbCorr (version 0.1.22; Küng P, 2023)
  • see (version 0.9.0; Lüdecke D et al., 2021)
  • tibble (version 3.2.1; Müller K, Wickham H, 2023)
  • R (version 4.4.1; R Core Team, 2024)
  • openxlsx (version 4.2.7.1; Schauberger P, Walker A, 2024)
  • ggplot2 (version 3.5.1; Wickham H, 2016)
  • forcats (version 1.0.0; Wickham H, 2023)
  • stringr (version 1.5.1; Wickham H, 2023)
  • rvest (version 1.0.4; Wickham H, 2024)
  • tidyverse (version 2.0.0; Wickham H et al., 2019)
  • readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  • dplyr (version 1.1.4; Wickham H et al., 2023)
  • purrr (version 1.0.2; Wickham H, Henry L, 2023)
  • readr (version 2.1.5; Wickham H et al., 2024)
  • xml2 (version 1.3.6; Wickham H et al., 2023)
  • tidyr (version 1.3.1; Wickham H et al., 2024)
  • knitr (version 1.49; Xie Y, 2024)
  • kableExtra (version 1.4.0; Zhu H, 2024)
report::cite_packages()